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We present a modified Lanczos algorithm to diagonalize lattice Hamiltonians with dramatically 
reduced memory requirements, without restricting to variational ansatzes. The lattice of size N 
is partitioned into two subclusters. At each iteration the Lanczos vector is projected into two 
sets of risvd smaller subcluster vectors using singular value decomposition. For low entanglement 
entropy 5*66, (satisfied by short range Hamiltonians), the truncation error is expected to vanish 
as exp(—ng^^ ''''). Convergence is tested for the Heisenberg model on Kagome clusters of 24, 30 
and 36 sites, with no lattice symmetries exploited, using less than 15GB of dynamical memory. 
Generalization of the Lanczos-SVD algorithm to multiple partitioning is discussed, and comparisons 
to other techniques are given. 

PACS numbers: 05.30.-d, 02.70.-c, 03.67.Mn, 05.50. +q 



I. INTRODUCTION 

Numerical ("exact") diagonalizations (ED) of quan- 
tum many-body Hamiltonians on finite clusters are of- 
ten used to advance our understanding of larger lattices. 
For example, Contractor Renormalization [1-3 uses ED 
to compute the short range interactions of the effective 
hamiltonian. ED on various size clusters [4] are indis- 
pensable as unbiased tests of mean field theories and 
variational wavefunctions. They are also used to obtain 
short wavelength dynamical correlations [5 and Chern 
numbers of Hall conductivity [6 . 

ED commonly use Lanczos algorithms [3 E], to effi- 
ciently converge to the low eigenstates. However, for a 
lattice of size N ^ with m states per site, the dimension 
of the Lanczos vectors (which are stored in the dynam- 
ical memory) increases as . Therefore, ED on larger 
lattices are prevented primarily by memory limitations, 
rather than processor speed. 

The central idea of this paper is to significantly reduce 
the memory cost, in order to enable ED of larger lattice 
sizes. We use singular value decomposition (SVD) to 
compress all Lanczos vectors into sets of 2nsvd vectors of 
size m^/^. 

As long as entanglement entropy of the target eigen- 
states obeys See « N/2\og{m) [9 , one can greatly 
economize on memory while maintaining high numerical 
accuracy. Many of the important many body Hamilto- 
nians of condensed matter (e.g. Hubbard and Heisen- 
berg models) have short range interactions. As a con- 
sequence, their ground states possess low entanglement 



*This work was supported by the U. S. DOE, Contract No. DE- 
AC02-76SF00515. 



entropy [TQHT3]. 

The idea of exploiting low entanglement entropy to 
compress wavefunctions by SVD, is not new. This is the 
key to the remarkable success of density matrix renor- 
malization group (DMRG) which has been used ex- 
tensively to obtain low energy state correlations of a 
large variety of Hamiltonians. DMRG is equivalent to 
variational minimization in the space of matrix prod- 
uct states [TOl [151 [IS^ • Extensions to wavefunctions with 
longer range correlations were given by multiscale entan- 
glement renormalization [TT. Nevertheless, sequential 
minimization may sometimes get "stuck" in false min- 
ima and not converge to the ground state. Therefore 
Lanczos methods are often called for to independently 
test the variational results. 

The paper is organized as follows. We begin by defin- 
ing the SVD for a bipartite split of the lattice, and pro- 
ceed to explain in practice how perform a single Lanczos 
step followed by an SVD projection which prevents the 
expansion of the memory cost. We describe the interme- 
diate matrix manipulations needed for orthonormaliza- 
tions and diagonalizations. In Section |III[ we estimate 
the SVD truncation error after projection, as a function 
of nsvd- We relate the error estimation to the bipar- 
tite entanglement entropy Seei using a generic asymp- 
totic form for the entanglement spectrum which is based 
on a classical gas model. We test in detail, the conver- 
gence of the Lanczos-SVD algorithm for the spin half 
Heisenberg antiferromagnet on Kagome clusters of up to 
36 sites. The entanglement spectrum asymptotics are 
verified for a partitioning of a 30 site cluster. Fig. |3] The 
ground state energy of 36 sites with Lanczos-SVD con- 
verges to relative errors of 6 x 10~^ for the three star 
line. Fig. [2] and 1.3 x 10~^ for the three star triangle. 
Fig. |4j Here we use a desktop computer with less than 
15GB of memory and no lattice symmetries. The results 
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FIG. 1: Partitioning 24 sites Kagome cluster for application of FIG. 3: Partitioning of 30 sites containing 5 hexagons, 
the Lanczos-SVD algorithm. The ground state entanglement Ground state entanglement entropy is See ~ 1.27. 
entropy is See — 1.51. 
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FIG. 4: Partitioning 36 sites. Ground state entanglement 
entropy is See ~ 2.59. 



II. THE LANCZOS-SVD STEP 

FIG. 2: Partitioning 36 sites for three stars in a line. Ground 

state entanglement entropy is « 1.12. ^_ Bipartite Singular Value Decompositi. 



agree with our analytical estiraate of the truncation error. 
In Section [VT| we discuss a possible extension of this ap- 
proach to multi-partitioning, and estimate the optimal 
reduction in memory cost that could be achieved. We 
conclude by a discussion which elaborates on the relative 
advantages and disadvantages of Lanczos-SVD, standard 
Lanczos and DMRG. 



Any state |?/^) of the full cluster (see as an example 
Fig. [T]) can be represented in a unique SVD form as 

OL 

= 1, (a|a')i = <5„a', (1) 

OL 

where the Aq, are positive, and \oi)i are "small" basis 
vectors of subclusters z = 1,2 (the subclusters on the 
two sides of the partitioning). Truncating the sum into 
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the largest n^vd terms defines the SVD projector, 



■ svd I 



(2) 



which introduces a wavefunction error e = Y^^^^ A?,. 



B. Application of PsvdH 

Lanczos-SVD economizes on the storage space by ap- 
plying an SVD projection after each application of the 
Hamiltonian on the Lanczos vector, 



PsvdH[ 



(3) 



The projection entails the following computational steps. 
H can be written as a sum of products of the two sub- 
cluster operators. 



M 



(4) 



where includes all internal interactions of subcluster 
i. (g) H2 is a product of operators residing on both 
subclusters. For example, a nearest neighbor Heisenberg 
model C^ij Si • Sj) with K bonds connecting the two 
subclusters has M = 2 + 3K terms. For example in 
Fig. [4l = 7 and M = 23. 

Acting with H on |?/^) produces a new state. 



HI 



nsvdM 

E 



(5) 



where the new (non orthonormal) small vectors are la- 
beled by the fused index u = (/i, a) , i.e. |z/)i = H^\a)i. 
The state ([5| lies outside the n^vd subspace, and we need 
to project it back using Pgvd in order not to further ex- 
pand the memory cost. We first ortho normalize by 
diagonalizing the Hermitian overlap matrices (of row di- 
mensions TlsvdM) 



1/3). 



E 



i = 1,2 



(6) 



where Di are diagonal and positive semidefinite, and Vi 
are unitary. = 1,2 are orthonormal sets in their 

respective spaces. Thus, the new vector is given by 



E 



Now we perform an SVD on the matrix C, 
C = plUlk'U2, 



(7) 



(8) 



where Vc is the normalization. Ul^Ul are unitary matri- 
ces which diagonalize the Hermitian products CC^ and 
C^C respectively. After computing A, /7i,/72 we obtain 
the SVD form of the new state. A is diagonal and nor- 
malized to TrA^ = 1 with positive eigenvalues A'^. We 
keep only the ngvd largest A'^, and obtain 



P,vdif|^) = E^"l")'ll")2' 



a=l 

where the new small vectors of i = 1, 2 are. 



^ — ■ V / au 



(9) 



(10) 



III. ENTANGLEMENT SPECTRUM AND SVD 
TRUNCATION ERROR 

The Lanczos algorithm rotates a set of basis states 
\ipn) into the lowest energy eigenstates with which the 
basis has a finite overlap. If we choose e to be much 
smaller than the lowest relative energy gap, the Lanczos- 
SVD vectors converge to a states which are within e dis- 
tance from the SVD projection of the corresponding ex- 
act eigenstate. 

To get an idea of how e(nsvd) converges, we must know 
the "entanglement spectrum" {s^} defined by A^ = e~*". 
Sa are pseudo-energies of the entanglement spectrum. A 
generic density of states can be modeled by a power law 
form 



r(p + i)' 



p> -1, 



(11) 



which describes the many-body density of states of a clas- 
sical gas with constant (Dulong- Petit) specific heat [18]. 
p counts with the number of entangled degrees of free- 
dom. The corresponding entanglement entropy is easy to 
evaluate, 



5'ee = - V A^ log(A^) = / dsspp{s)e-' = p + 1. 
Choosing a high cut-off exponent Sc such that 



(12) 



^svd 

we arrive at the error estimate 
dspp{s)e~ 



r(p + i) 



(13) 



(14) 



Combining Eqs. (13) and (14), yields for See >> 1, the 



asymptotic expression 



(15) 



Hence by choosing the ratio n^vd/^^'"'' » I one ensures 
an exponentially small truncation error. 
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FIG. 5: Entanglement spectrum. We extract the asymptotic 
behavior of the entanglement density of states pp, for the first 
excited state of the five hexagon cluster of Fig. [3] The line 
depicts a pure power law at large pseudoenergies s, consistent 
with the Dulong-Petit form given in Eq. (11). 
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FIG. 6: Lanczos-SVD truncation errors, for the 24 sites 
Kagome cluster of Fig. ^ risvd is the number of retained 
SVD states of the ground state ([3|. e is the wavefunction 
error Eq. ^ and following text). AEq/Eq is the relative er- 
ror in the the Lanczos-SVD ground state energy as compared 
to the exact (standard Lanczos) result. The rapid decay of 
the errors for low values of risvd << 2^^ is due to the low 
entanglement entropy, see Section [ill] 



basis" of dimension L 



n+L 

E 



(16) 



In Figure |5j we plot the entanglement spectrum for 
the first excited eigenstate of the 30 site, five hexagon 
Kagome cluster depicted in Figure |3] The log-log plot 
demonstrates the asymptotic power law density of states 
~ 5^, where p ~ 2.55. We note that for this system 
we are in the low See regime and hence the difference 
between the fitted value of p + 1 and the entanglement 
entropy of 1.27. Nevertheless, the density of states at 
high pseudo-energies extrapolates well to the asymptotic 
power law behavior of a classical gas as modelled in this 
section. 



IV. IMPLEMENTATION OF THE 
LANCZOS-SVD ITERATION 



where A are the coefficients determined by diagonalizing 
the overlap matrix (see e.g. Eq. (|6|). 

Subsequently, we compute the matrix elements of the 
reduced Hamiltonian, 



1,...L. 



(17) 



The reduced Hamiltonian matrix is diagonalized. Its 
lowest eigenvalue and eigenvector yield the best approxi- 
mation to the ground state at this level of iteration [19] . 
We bring the resulting wavefunction to the SVD form, 
again truncated into ngvd terms. It becomes the new ini- 
tial state I '0(^+^+1)) for the next Lanczos iteration. Ex- 
cited states can be calculated by starting with an initial 
state orthogonal to the converged lower energy states. 



The Lanczos-SVD routine proceeds as follows: We ini- 
tialize \ip)^^'> as a direct product of the two subcluster 
states. We compute (T^svci^j^'IV^^^^) = IV^^""^) as de- 
scribed above. Since our method is economical in mem- 
ory, we can afford to retain L sequential Lanczos vectors 
I'^C^))^ 1^0(^+1))^ . . . 1?/;'^^+^)), which speeds up the conver- 
gence with iteration number considerably. (If memory is 
scarce, one could use the slower method of keeping only 
two Lanczos vectors). 

Now, we compute the overlap matrix and orthonormal- 
ize this set of Lanczos vectors. This produces a "rotating 



V. NUMERICAL TESTS 

For the spin-half Heisenberg antiferromagnet 

H = Y,S,-S,, (18) 

(ij) 

we tested the convergence of the Lanczos-SVD algorithm 
on Kagome clusters depicted in Figs, [l] [2j [3]and[4j 

In Fig.[6j we plot the Lanczos-SVD truncation error for 
the ground state wave function and energy, as a function 
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P = 8, p = 3 



-5 



bJO " 
O 



-9 



\ 3 Star Triangle 












\ ° 


X^Star Line 




\ -10 






\ -20 






) 50 100 150 


20 


40 


60 



iteration i 




FIG. 8: Multiple subclusters 



FIG. 7: Lanczos-SVD error convergence as a function of iter- 
ation, for 36 site Kagome clusters. The main plot is for Fig.[4j 
and the inset is for Fig. |2] AEo{i) is difference between the 
energy at iteration i and the converged ED result given by 
standard Lanczos [20^. For both clusters we used risvd — 200 
and L = 4 Lanczos vectors. 



of nsvd- The errors decrease rapidly on the logarithmic 
scale as expected by Eq. (??), arriving at ~ 10~^^ for 
nsvd = 200. 

For the 30 site Kagome cluster depicted in Fig ??, the 
energies of the four lowest S=0 eigenstates, and the first 
triplet S=l state, converged to a very high accuracy of 
10-^^ using nsvd = 200. 

In Fig. [7| we show the convergence of ground state 
energy of Lanczos-SVD versus iteration for clusters of 36 
sites. We use nsvd = 200, and L = 4. For the three star 
triangle, the exact ground state energy as determined 
by standard Lanczos is £^o(36sites) = —14.859397 [20 . 
The entanglement entropy is See ^2.5. The calculation 
converges to relative energy accuracy of 1.3 x 10~^, and 
an SVD truncation error of similar magnitude. 

In the inset of Fig. [7| we show a much smaller error 
for the three star line of Fig. [2j which converged to a 
relative error of 6.3 x 10~^. This is to be expected since 
the entanglement entropy of the linear arrangement of 
the three stars is only See ^ 1.12. 

The numerical tests were therefore consistent with Eq 



(15). 



All the above calculations were performed using multi- 
core workstations. The maximum memory usage was 
kept under 15 GB of memory, even though no lattice 
symmetries were implemented in the computations. The 
time required for the most intensive calculations (36 sites, 
nsvd = 200) on using parallelisation with 16 cores was a 
little more than 70 minutes per iteration Fig |4] and about 
35 minutess for Fig. [2] A serial MAPLE 15 implemen- 
tation used in the computations for the 30 site case took 
about one day for each eigenstate for the same nsvd- 



VI. EXTENSION TO MULTIPLE 
PARTITIONING OF LARGE LATTICES. 

The Lanczos-SVD compresses the memory require- 
ment by a single division of the cluster into two subclus- 
ters i = 1,2. This idea could be extended to recursive 
partitioning [17]. For the sake of crude memory estima- 
tion, each small vector (e.g. \a)i in Eq. (TJ can be de- 
composed into nsvd products of even smaller subcluster 
vectors. If the SVD is thus iterated p times, one obtains 
a representation in terms of small vectors of P = 2^ sub- 
clusters, is thus stored in terms of a set of the smallest 
vectors. The memory cost after applying a Lanczos step 
is as follows. 

For concreteness, let us consider a two dimensional disk 
of radius >> 1, containing TV ~ ttR'^ sites of spin half, 
divided into P equal sections as shown in Figj8| The 
sections are labeled by a binary number i = (^i, ^2, • • • v)' 
ik = 0, 1. The recursive SVD decomposition yields the 
expression 



E 

a.1 ,q;2 . . .( 



n 



,ap)i. 



(19) 

The SVD weights are labeled according to the bound- 
aries they describe, as shown in Fig. ([s]). Each ai runs 
over nsvd numbers, which means that each section i is 
represented by n^^^ vectors of dimension 2^/^. By the 
"area law" See oc on each boundary. As shown before, 
we must retain nsvd ^ terms in each SVD, where 
c(e) > 1 in order to achieve a desired truncation error e. 

After applying to we generate a factor of M ^ 
6R more small vectors. Thus we should store 6PRn^^^ 
small vectors. Thus the memory cost is 

7tR^ log(2)\ 



Mc ~ 6Pi?exp (cR\og{P) 



(20) 



Minimizing A4c{P) one finds the optimal partitioning 
P^P^^ and the optimal memory cost A4^p^ at large N to 
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scale as 

popt _ 7rl0g(2)J? ^ (\/]vAlog(Af/7r)-2) 

(21) 

This would amount to a significant compression of mem- 
ory as compared to standard Lanczos A4 ^ 2^ . The re- 
maining challenge is to speed up the significantly larger 
computational time needed to orthonormalize and SVD 
large sets of small vectors. 

VII. DISCUSSION 

In this paper we have exploited the SVD to expand 
the lattice sizes which can be treated by the Lanczos 
algorithm for ED. Relative to standard Lanczos, Lanczos- 
SVD demands longer computation time due to additional 
matrix manipulations. 

Does the SVD projection interfere with the Lanczos 
convergence?. The SVD projection introduces a trun- 
cation error in the rotating vector. However, the Lanc- 
zos vector is rotated toward the ground state (or some 
other target state) when the relative energy splitting to 
neighboring eigenstates is larger than the truncation er- 
ror. This rotation ceases once the energy has converged 
within the SVD error. Increasing ngvd will allow further 
convergence. It is simplest to think about the SVD pro- 
jection in the same footing as an additional fioating point 
error which limits the numerical convergence of standard 
Lanczos algorithms. In computing Fig. [7|we have indeed 
verified that the implementation of SVD projection in 
each iteration, does not slow down the energy conver- 
gence per each Lanczos step. 

Computing time of Lanczos-SVD may be significantly 
reduced by using multiple cores and by parallelizing the 
code. In particular the overlap calculations, which are 
currently the most time-consuming, are easily parallelize- 
able. In standard Lanczos, memory reduction can be 
achieved by exploiting lattice and spin symmetries, which 
demands special boundary conditions and extensive pro- 



gramming. Lanczos-SVD is therefore simpler to imple- 
ment, and can address arbitrary boundary conditions. 

It is often asked what advantage a Lanczos based ED 
has over DMRG and related variational methods? The 
answer of course depends on the purpose of the calcu- 
lation. DMRG has proven very efficient (especially for 
ground state correlations) for larger lattices than Lanczos 
methods can address. However, in practical applications, 
DMRG advances toward the ground state by sequential 
minimizations, e.g. "sweeping" the parameters of the 
wave function in real space. In cases of high frustra- 
tion and competing phases (e.g. near a phase transition), 
sweeping methods can get stuck in metastable states [21] . 
(Consider e.g. the difficulty of getting rid of defects in a 
phase separated system by sweeping methods). 

The Lanczos step, on the other hand, does not neces- 
sarily move the state in the direction of maximal slope. 
Instead, if numerical accuracy is sufficient, it steadily ro- 
tates the Lanczos vector toward the true ground state (or 
some other target eigenstate). 

Lanczos-SVD can therefore be used to determine the 
ground state and low excitations of limited size clusters 
with well controlled accuracy. It could be used to check 
DMRG convergence and test variational ansatzes. As 
mentioned in the introduction, a primary purpose for us- 
ing ED on small clusters is to derive an effective Hamil- 
tonian by the CORE method [22 . The CORE effective 
Hamiltonian can then be studied on the coarse grained 
lattice by iterating CORE, or by variational methods. 
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